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The widespread use of antibiotics in association with high-density clinical care has driven the emergence of drug- 
resistant bacteria that are adapted to thrive in hospitalized patients. Of particular concern are globally disseminated 
methicillin-resistant Staphylococcus aureus (MRSA) clones that cause outbreaks and epidemics associated with health care. 
The most rapidly spreading and tenacious health-care-associated clone in Europe currently is EMRSA-15, which was first 
detected in the UK in the early 1990s and subsequently spread throughout Europe and beyond. Using phylogenomic 
methods to analyze the genome sequences for 193 S. aureus isolates, we were able to show that the current pandemic 
population of EMRSA-15 descends from a health-care-associated MRSA epidemic that spread throughout England in the 
1980s, which had itself previously emerged from a primarily community-associated methicillin-sensitive population. The 
emergence of fluoroquinolone resistance in this EMRSA-15 subclone in the English Midlands during the mid-1980s 
appears to have played a key role in triggering pandemic spread, and occurred shortly after the first clinical trials of this 
drug. Genome-based coalescence analysis estimated that the population of this subclone over the last 20 yr has grown 
four times faster than its progenitor. Using comparative genomic analysis we identified the molecular genetic basis of 
99.8% of the antimicrobial resistance phenotypes of the isolates, highlighting the potential of pathogen genome se¬ 
quencing as a diagnostic tool. We document the genetic changes associated with adaptation to the hospital environment 
and with increasing drug resistance over time, and how MRSA evolution likely has been influenced by country-specific 
drug use regimens. 

[Supplemental material is available for this article.] 


Methicillin-resistant Staphylococcus aureus (MRSA) remains a ma¬ 
jor global cause of health-care-associated infections since it was 
first described five decades ago following the introduction of 
penicillinase-stable (3-lactam antibiotics into clinical practice 
(Klein et al. 2007). During this period, MRSA strains have emerged 
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several times by acquiring variants of staphylococcal cassette 
chromosome mec (SCC mec) elements that carry the mecA methi- 
cillin resistance determinant (Robinson and Enright 2003; Niibel 
et al. 2008). However, the vast majority of MRSA isolated world¬ 
wide belong to a limited number of clones, some of which are as¬ 
sociated with global epidemics. 

Epidemic MRSA-15 (EMRSA-15) has proven to be particularly 
successful over the last two decades, transmitting rapidly within 
and between hospitals, as well as to different countries. EMRSA- 
15 carries a type IV SCC mec element and belongs to multilocus 
sequence type ST22. Initially isolated in the southeast of England 
in 1991 (Richardson and Reith 1993), the spread of EMRSA-15 in the 
UK health-care setting was rapid; by 2000, EMRSA-15 accounted 
for over 60% of MRSA nosocomial bacteremias in England 
(Johnson et al. 2001). Over this period the proportion of MRSA 
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among S. aureus bacteremia in the UK increased from <2% to 40% 
(Johnson et al. 2001). 

Shortly after EMRSA-15 was reported in England, MRSA in¬ 
distinguishable from EMRSA-15 were reported to cause outbreaks 
in other geographic regions, including Germany, the Czech Re¬ 
public, Portugal, New Zealand, Australia, and Singapore (Pearman 
et al. 2001; Witte et al. 2001; Hsu et al. 2005; Melter et al. 2006; 
Amorim et al. 2007; Grundmann et al. 2010). These belong to 
multilocus sequence type ST22 or its close relatives. In a recent 
structured survey, ST22 was the most frequently found MRSA 
across Europe (Grundmann et al. 2010), and it was among the 
three most frequent MRSA clones in seven out of 15 countries for 
which quantitative data were available (Grundmann et al. 2010). 
Although reported in many countries, ST22-MRSA accounted 
for only 0.2% of MRSA in a surveillance study recently performed 
in the USA (Limbago et al. 2009) and has not, as yet, been re¬ 
ported from South America (Sola et al. 2008) and from large 
parts of Asia (Arakere et al. 2005; Ko et al. 2005). ST22-MRSA has 
repeatedly demonstrated its ability to supplant and replace other 
formerly established MRSA strains (Hsu et al. 2005; Amorim 
et al. 2007; Aires-de-Sousa et al. 2008; Witte et al. 2008). Of fur¬ 
ther concern, ST22-MRSA containing the lukS/F genes (encoding 
the Panton-Valentine leucocidin) has been reported to cause 
community-associated infections (Linde et al. 2005; Witte et al. 
2007). 

The mechanisms that lead to the emergence and spread of 
MRSA strains are of great interest. Despite intense research in this 


area, the spatiotemporal dynamics of the population structure and 
epidemic spread remain poorly understood, partially due to the 
limited resolution provided by contemporary genotyping tools. 
In this study we have investigated the origins and evolution of 
EMRSA-15 by whole-genome sequencing 193 ST22 isolates from 
15 countries, collected between 1990 and 2009. Integrating evo¬ 
lutionary and spatial information, we have reconstructed the spa¬ 
tial and temporal dynamics underpinning the expansion of this 
clone and ascertained the genetic changes correlating with its en¬ 
hanced spreading success. 

Results 

Phylogenomic reconstruction of the ST22 lineage 

We determined full-genome sequences from 193 ST22 isolates 
(Supplemental Table SI) including isolates both from hospital and 
community settings, and also 12 isolates from the earliest de¬ 
scription of EMRSA-15 (Richardson and Reith 1993). Once mobile 
genetic elements (MGEs) were excluded, a total of 8095 single 
nucleotide polymorphisms (SNPs) in the "core" genome were used 
to reconstruct the evolutionary history of S. aureus ST22 (Fig. 1). 
These data revealed that the core genome has accumulated varia¬ 
tion at a rate of 1.3 X 10 -6 substitutions per nucleotide site per year 
(95% Bayesian credible intervals, 1.2 X 10 -6 to 1.4 X 10 -6 ; Sup¬ 
plemental Table S2). Substitution rates varied slightly among dif¬ 
ferent phylogenetic clades within ST22, but these differences 
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Figure 1. Phylogeny of ST22 and the emergence of MRSA clones. ( A ) Maximum likelihood phylogenetic tree of ST22 isolates. The tree was rooted by 
using the distantly related 5. aureus isolate MSSA476 as an outgroup. Colors indicate the isolates' countries of origin. Roman numerals indicate acquisitions 
of structurally different SCCmec elements, which cause methicillin resistance. ( B ) Maximum clade credibility tree of the ST22-A clade based on BEAST 
analysis using a variable clock rate (uncorrelated lognormal) model. Tips of the tree are constrained by isolation dates, the time scale is shown at the 
bottom. Gains and losses (A) of genetic determinants for resistance to methicillin (SCCmec IVh), fluoroquinolones (point mutations in grIA and gyrA), 
erythromycin (plasmid-encoded ermC), and clindamycin (mutations in ermC leader peptide region, c-ermC ) have been mapped on the tree by applying 
the parsimony criterion. 
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were not significant (Supplemental Table S2). Our rate estimate 
was similar to those previously reported for other Staphylococcus 
aureus lineages (Lowder et al. 2009; Harris et al. 2010; Niibel 
et al. 2010). This estimate then provided a temporal calibration 
with which to date the emergence and geographic spread of 
EMRSA-15, and evolutionary events that have contributed to its 
success. 

The level of homoplasy in the SNP data set was extremely 
low (homoplasy index, 0.009), reflecting the clonal propagation 
of the core genome. The majority of isolates (162 of 193, 84%), 
collected from multiple countries, clustered in a single clade 
(ST22-A) (Fig. 1A), which encompasses the epidemic strain EMRSA-15. 
All but four of the ST22-A isolates carried a single nucleotide de¬ 
letion within the ureC gene, resulting in a frame-shift inactivation 
of the urease alpha subunit and the consequent nonproduction 
of urease (Supplemental Table SI), which is a phenotypic char¬ 
acteristic of EMRSA-15 (Pearman et al. 2001). The four isolates 
lacking this deletion are the most basal isolates of the ST22-A 
clade (for a list of the genetic variation that distinguishes the 
ST22-A clade from the rest of the ST22 population, see Supple¬ 
mental Table S3). 

The origin of 173 of the isolates was recorded. Whereas 
the vast majority (97%) of ST22-A isolates were collected in 
hospitals, this proportion was significantly smaller in non- 
ST22 isolates (P < 0.001; Supplemental Table SI). This sug¬ 
gests that members of ST22-A have adapted to transmission 
and survival in health-care settings. In contrast, other clades 
within ST22 encompass methicillin-susceptible S. aureus iso¬ 
lates, MRSA with different types and subtypes of SCC mec ele¬ 
ments, and isolates with genes for the Panton-Valentine leuco- 
cidin toxin. The phylogeny indicates that these isolates have 
distinct evolutionary origins from the health-care-associated 
clone ST22-A isolates (Supplemental Table SI) and represent the 
ST22 founding population from which ST22-A/EMRSA-15 re¬ 
cently emerged. 

Included in the sequenced collection are 12 isolates from the 
original description of EMRSA-15 (Richardson and Reith 1993). 
At the time, Richardson and Reith noted that the majority of 
EMRSA-15 in the UK were ciprofloxacin sensitive; however, they 
identified a group of variant isolates from the English West Midlands 
(Birmingham region) that were ciprofloxacin resistant. All eight 
ciprofloxacin-sensitive EMRSA-15 isolates from the Richardson 
and Reith study are found in the ST22-A1 clade (Fig. IB; Supple¬ 
mental Table SI), whereas the four ciprofloxacin-resistant isolates 
form a basal group in clade ST22-A2. It therefore appears that ST22- 
Al, while spreading among hospitals in the British Isles, formed 
a progenitor population from which a fluoroquinolone-resistant 
pandemic clone (i.e., ST22-A2) emerged. 

EMRSA-15 genotype 

A combination of gene acquisition and core genome muta¬ 
tion shapes the genome of ST22. Therefore, in addition to SNP 
discovery based on mapping sequencing reads to a reference 
genome, we also determined the full genetic repertoire from 
each of the isolates, including MGEs, by de novo assemblies of 
sequencing reads. Using the phylogenetic reconstruction, we 
have been able to examine the evolution of ST22-A, and also 
pinpoint the genetic events that accompanied the emergence of 
EMRSA-15. 

A significant component of the S. aureus genome is derived 
from MGEs that contribute to the accessory genome and includes 


SCC mec, prophages, S. aureus Pathogenicity Islands (SaPI), IS ele¬ 
ments, transposons, and plasmids (Lindsay and Holden 2004). 
Comparative analysis reveals the composition and fluidity of the 
ST22 accessory genome (Fig. 2), and this accounts for between 4% 
and 14% of CDS within the ST22 genomes. Comparison of the 
ST22-A and non-ST22-A strains reveals that, typically, the accessory 
genome of ST22-A strains is larger and less variable than that of non- 
ST22-A strains (Fig. 2A). 

Analysis of the origins and location of CDSs in the ST22 acces¬ 
sory genomes show MGEs belonging to a range of classes (Fig. 2B), 
with prophages being the largest component. The accessory ge¬ 
nome variation within the ST22-A isolates results from both the 
loss of ancestrally acquired elements such as prophages and SCC mec 
(Fig. 2) and the gain on multiple occasions of plasmids and pro¬ 
phages (Supplemental Table SI). 

The ST22 reference genome, HO 5096 0412, is a representa¬ 
tive of EMRSA-15 and is found in the ST22-A2 clade (Supple¬ 
mental Table SI). Among the MGEs it contains are the following: 
an SCC mec IVh element (Milheirico et al. 2007), which carries the 
mecA gene encoding methicillin resistance; a Tn554 like-transposon, 
encoding a (3-lactamase; two prophages (4>Sa2 H o 5096 0412 / which 
encodes a Bcgl-like restriction enzyme, and cj>Sa3Ho 5096 0412 , 
which encodes staphylokinase [Sako and Tsuchida 1983]) and 
the immune evasion proteins CHIPS and SCIN (van Wamel et al. 
2006); and a 2.4-kb ermC carrying plasmid (Koser et al. 2012). 
Comparative genomic analysis of the ST22 isolates showed that 
these elements were conserved in the ST22-A population, with 
few exceptions (Supplemental Table SI). 

Using the phylogeny we showed that the SCC mec IVh ele¬ 
ment was acquired once, just prior to the emergence of ST22-A. 
Similarly, c|>Sa2 appears to have been acquired prior to the 
emergence of ST22-A (Supplemental Table S3). Prophage 4>Sa3 
seems to have been acquired before the emergence of the ST22-A 
lineage because it is present in members of the non-ST22-A 
population. A SaPI element encoding enterotoxin C (Novick 
2003) is absent from the reference genome, but present in —68% 
of the ST22-A isolates (Supplemental Table SI). This element is 
also found in non-ST22-A isolates at the same site (3' region of 
GMP synthase, SAEMRSA1503430), suggesting that like c|>Sa3, 
its association with the ST22 lineage predates the emergence 
of ST22-A. In addition, other MGEs are variably present in the 
pandemic clone, several of which are associated with the 
transfer of antibiotic resistance genes; but, unlike the elements 
above, none of these have been stably maintained since ST22-A2 
emerged. 

In addition to the MGE-induced changes in the genomic 
architecture of ST22-A, all ST22-A genomes possess five indels in 
the core genome, including a 2268-bp deletion in the fibronectin- 
binding protein (FnBP) locus (Supplemental Table S3). This FnBP 
deletion was caused by homologous recombination between the 
C-terminal region of fnbA and the N-terminal region of fnbB, 
resulting in a gene fusion, fnbA-fnbB. FnBPA and FnBPB bind the 
host extracellular-matrix proteins fibronectin, fibrinogen, and 
elastin (Greene et al. 1995; Roche et al. 2004), allowing the 
S. aureus to connect to cellular integrins and trigger internaliza¬ 
tion (Fowler et al. 2000). In the ST22 isolates containing fnbA and 
fnbB (non-ST22-A), these two homologs are —71% identical 
(—67% identical at the amino acid level) and 2871 bp and 3045 bp 
in length. The predicted recombination site occurs —1850 bp into 
fnbA, thereby generating a gene fusion that retains most of the 
sequence of fnbA. The fusion protein is the same length as FnBPA 
and is 98.5% identical, whereas it is only 68.6% identical to 
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Figure 2. Accessory genome diversity in ST22. CDSs in the accessory genome of each isolate were clustered and assigned to homology groups (HGs). 
(A) Conservation of accessory genome HGs across the ST22 population. The isolates are ordered according to the phylogeny displayed along the left and 
top of the figure. The conserved HGs in the pairwise comparison are displayed as a heat map matrix and colored according to the number of conserved HGs 
(range 64 to 235 HGs) per pair; see figure for scale. ( B ) Size and composition of each isolate's accessory genome based on the number of HG CDSs. The HGs 
have been subdivided into different mobile genetic element types based on matches to an annotated MGE reference set included in the clustering (see 
legend to figure). 


FnBPB. It is therefore likely that the encoded fusion protein 
retains most of the binding properties of FnBPA (Roche et al. 
2004). 

From the phylogeny we defined the genetic changes that 
differentiate the ST22-A2 and ST22-A1 populations (Supplemental 
Table S4). That differentiation was not marked by any unique 
acquisitions of MGEs. Instead, the variation in the core genome 
that marks this split includes a single insertion and 30 SNPs, of 
which 18 result in nonsynonymous substitutions, and one non¬ 
sense mutation. The nonsense mutation introduces a premature 
stop codon in the ECM-binding protein homolog Ebh. This very 
large surface protein binds fibronectin and contains a C-terminal 
transmembrane domain that anchors the intact protein in the 
cell envelope (Clarke et al. 2002). The truncated variant caused 
by the nonsense mutation lacks this transmembrane domain, 
and ST22-A2 isolates likely do not express this adhesin on their 
surface. 


Demography and phylogeography 

The high proportion of ST22-A isolates in the sample and their 
very low sequence diversity, points to ST22-A being a highly fit, 
recently emerged MRSA strain. This inference is further sup¬ 
ported by coalescence analyses estimating that the popula¬ 
tion size of ST22-A has grown exponentially by 24% each year 
(Bayesian credibility intervals, 20%-29%) since the 1990s, which 
is four times faster than the average growth of other ST22 strains 
(non-ST22-A) (Fig. 3). The rapid population expansion is accom¬ 
panied by the star-shaped branching pattern of the ST22-A clade 
(Fig. 1A). 

A Bayesian phylogeographic reconstruction of ancestral 
nodes in the phylogenetic tree supports a root of the ST22-A clade 
in the UK. Within a few years of its emergence there, multiple ex¬ 
ports spread it from the UK to Australia as well as to various coun¬ 
tries in Europe, Asia, and Africa (Fig. IB; Supplemental Fig. SI). 
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Figure 3. Exponential growth rates (X) of ST22-A and other ST22 
isolates (posterior probability density functions from Bayesian analy¬ 
sis). ST22-A is estimated to have grown significantly faster than 
non-ST22-A. The population size N at time t was modeled by N(t) = 
N( 0)xe (xx °, where X is the exponential growth rate and t is measured 
in years. Accordingly, the population size change per year is e\ 


Early ST22-A isolates (clade ST22-A1 in Fig. IB) were susceptible to 
fluoroquinolones, and apparently did not spread beyond the 
British Isles. Bayesian modeling of the spread of fluoroquinolone- 
resistant ST22-A2 in the UK suggests that ST22-A2 originated 
within the English Midlands (Fig. 4), and was restricted to this 
geographical region until the late 1980s. The same analysis further 
suggests that isolates from the Midlands seeded infection to the 
north, reaching Scotland by around 1993, but also to the south, 
reaching London around 1998 (Fig. 4; Supplemental Movie SI). By 
2000, the epidemic had spread to most locations within the UK 
via multiple routes and sources (e.g., Bristol, from the Southern 
England and Midlands; Manchester, from the Midlands and the 
North; Aberdeen from elsewhere in Scotland and North England). 
The UK epidemic then spread globally through multiple inter¬ 
national transmission events (Fig. 1). 

Estimates for the timings of first appearances of ST22-A in 
countries beyond the UK (Supplemental Fig. SI) are concordant 
with available epidemiological evidence for EMRSA-15, suggesting 
epidemic expansion in those countries soon after its introduction 
(Richardson and Reith 1993; Pearman et al. 2001; Witte et al. 2001; 
Hsu et al. 2005; Melter et al. 2006). Our data reveal multiple intro¬ 
ductions from the UK into Germany, Denmark, and Australia, pos¬ 
sibly via the migration of colonized health-care workers (Pearman 
et al. 2001). Country-specific clades in the phylogenetic tree point 
to subsequent diversification within many recipient countries, 
and the data support additional spreading events, e.g., from 
Portugal to Germany (Fig. IB). Because basal isolates outside the 
ST22-A clade are also geographically widespread, ST22 was almost 
certainly present in many countries prior to the emergence of 
ST22-A (Fig. 1A; Supplemental Table SI). Indeed, methicillin- 
susceptible ST22 is known to form a common component of the 
nasal carriage population in many places (Feil et al. 2003; Holtfreter 
et al. 2007; Ruimy et al. 2009). 


The development of drug resistance in ST22 

The broad range of resistance traits within ST22-A contrasts strik¬ 
ingly with the more susceptible and diverse population from 
which it emerged (SupplementalTable SI). For ST22-A, the number 
of resistance traits per isolate increased with time (P < 0.001) (Fig. 5), 
whereas the trend in the rest of the population (non-ST22-A) was 
not statistically significant (P > 0.1). 

In order to investigate the evolutionary events responsible 
for the increase in resistance, we analyzed the whole-genome se¬ 
quences to elucidate the antibiotic resistance genotype of ST22. 
Molecular determinants within the genomes explained 99.8% of 
the measured phenotypic resistance traits (847 resistance traits in 
total). Only two phenotypes did not correlate with the predicted 
genotypes. HO 7258 0475 05 was erythromycin resistant, but 
lacked ermC and all other erythromycin resistance mechanisms 
that have been documented in the literature. Strain 09-00678 con¬ 
tained thyA and dfrA, but was not phenotypically cotrimoxazole 
resistant. For both isolates it is possible that the drug-resistance 
determinants were lost during lab culturing prior to DNA prepara¬ 
tion or phenotypic testing. 

Our genomic analyses showed that evolution of resistance in 
ST22 has occurred in a stepwise fashion encompassing both the 
core and accessory genomes. Acquisition of a transposon-encoded 
(1-lactamase BlaZ confers penicillin resistance, and the widespread 



Figure 4. Bayesian reconstruction of the spread of ST22-A2 in the UK. A 
continuous spatial diffusion model was used to reconstruct the finer-scale 
geographical dispersal of ST22-A2 within the UK and to predict the origin 
of fluoroquinolone resistance. Lines indicate the inferred routes of spread 
with confidence displayed as green ovals representing 80% of the highest 
posterior density (HPD) for latitude and longitude. The timing of trans¬ 
mission events are represented by red (old) to black (recent) lines and 
light- to dark-green oval shading (for the animation of the reconstruction 
of the spread, see Supplemental File S2). The maps are based on satellite 
pictures made available at Google Earth (http://earth.google.com). 
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Figure 5. Increase in antibiotic resistance traits within ST22-A pop¬ 
ulation over time. Number of genotypic resistance determinants per 
ST22-A isolate (red, no. 162) and non-ST22-A isolates (blue, no. 31) from 
all countries, 1990-2008, including regression lines with 95% confidence 
interval for each group. Size of circle corresponds to number of isolates. 


distribution of this transposon in ST22 is consistent with a single, 
primordial development of (3-lactam resistance, and probably 
reflects usage of penicillin since the 1940s. Only six ST22 isolates in 
our collection are sensitive to penicillin, four of which have lost 
the Tn552-like transposon carrying the blaZ gene, and one of which 
carries a blaZ gene inactivated through an IS element insertion in its 
promoter region (Supplemental Table SI). The penicillin-resistant 
ST22 population then gave rise to several methicillin-resistant 
clones, which emerged through imports of SCC mec elements into 
their chromosomes. The phylogenomic distribution of different 
SCC mec variants in the phylogenomic context provides evidence 
of at least seven independent imports into ST22, one replacement 
of SCC mec type IV by SCC mec type V, and two secondary losses of 
SCC mec through precise excision by homologous recombination 
between flanking direct repeat sequences (Fig. 1; Supplemental 
Table SI). As described above, ST22-A emerged after the acquisition 
of SCC mec IVh, which is predicted to have occurred over a decade 
before the first report of EMRSA-15 (Fig. 1; Richardson and Reith 
1993). 

Around 1986 (95% Bayesian credible intervals 1984-1988) 
ST22-A developed fluoroquinolone resistance through two point 
mutations generating amino acid substitutions Ser80Phe in topo- 
isomerase IV (GrlA) and Ser84Leu in gyrase A (GyrA). This geno¬ 
type is found in the ST22-A2 clade (Fig. IB) that spread throughout 
the UK and globally. We note that our dating estimate for acqui¬ 
sition of fluoroquinolone resistance coincides closely with the in¬ 
troduction of fluoroquinolone drugs into clinical medicine in the 
UK, ciprofloxacin being first licensed in February 1987. Further, 
phylogeographic modeling suggested that fluoroquinolone resis¬ 
tance emerged in the English Midlands in the early 1990s (Fig. 4). 
Notably, clinical studies with ciprofloxacin occurred in the Mid¬ 
lands prior to the licensing of ciprofloxacin in 1987 (Crump et al. 
1983; Finch et al. 1986; Silverman et al. 1986). 

In the subsequent decades the pandemic has continued along 
a trajectory of increasing antibiotic resistance, developing re¬ 
sistance to a range of antibiotics of different classes (tetracy¬ 
clines, aminoglycosides, cotrimoxazole, fusidic acid, mupirocin, 
and oxazolidinones) on multiple occasions (Supplemental Fig. S2; 
Supplemental Table SI). The development of these resistances has 
been rather sporadic and also less stably maintained than those 


against (3-lactams and fluoroquinolones, but nevertheless dem¬ 
onstrates clear evidence of the strong selective pressure exerted 
by antibiotic usage and clinical practice. For example, our data 
set documents multiple independent events of erythromycin- 
resistance acquisition based on the uptake of plasmids encoding 
the ErmC erythromycin ribosomal methylase (Fig. IB). Phyloge¬ 
netic analysis also indicates that the ermC plasmid was lost mul¬ 
tiple times, and in several cases these losses affect entire clades, 
suggesting that they represent genuine events in the natural pop¬ 
ulation rather than laboratory artefacts (Fig. IB). ErmC addition¬ 
ally confers resistance to clindamycin when mutations in the 
leader peptide region of the ermC gene cause its constitutive ex¬ 
pression (Levin et al. 2005). All clindamycin-resistant isolates 
possess genetic rearrangements that disrupt the promoter region 
upstream of ermC, including deletions of part or all of the ermC 
leader peptide region or insertion of an IS element in ermC leader 
peptide region. At least eight independent, structurally different, 
mutations have generated clindamycin resistance in the pop¬ 
ulation examined (Supplemental Table SI) and five of these were 
found in ST22-A (denoted "c-erraC" in Fig. IB). The data also sug¬ 
gest that country-specific variation in antibiotic usage may be 
shaping the geographic distribution of resistance. The prescription 
level of clindamycin is much higher in Germany than in the UK 
(Coenen et al. 2006; Meyer et al. 2006; de With et al. 2011), causing 
selective pressure for a high prevalence of clindamycin resistance. 
Among the 34 erythromycin-resistant isolates from Germany, 32 
(94%) are also clindamycin resistant, whereas none of the 45 UK 
isolates are resistant to clindamycin (significant difference, P < 
0.001, Fisher's exact test). 

Discussion 

EMRSA-15 represents a "measurably evolving population" 
(Drummond et al. 2003), accumulating genetic variation over ep¬ 
idemiological timescales at an average rate similar to previous es¬ 
timates for other S. aureus strains (Lowder et al. 2009; Harris et al. 
2010; Nubel et al. 2010). One consequence of this fast nucleotide 
substitution rate is that evolutionary events could be dated, in¬ 
cluding the stepwise acquisition, refinement, and loss of multiple 
drug-resistance traits. Furthermore, past temporal changes in pop¬ 
ulation size and geographic expansion could be analyzed retro¬ 
spectively based on the molecular traces that these processes left in 
the genomes. These conclusions were possible because our genome 
sequences span a wide range of time and space, which is a pre¬ 
requisite for such coalescence-based analyses. Our data revealed 
much more detail about the pathogen's population history than 
would have been possible with a snapshot survey from a single 
time point. 

In principle, changes in the frequency of a genetic variant 
maybe due to stochastic events or genetic adaptation and selection 
(Keim and Wagner 2009). For the highly successful MRSA strains 
described here, it seems plausible that they are more specialized to 
the restricted ecological niches in the hospital environment, where 
exposure to antibiotics is common and survival depends on effi¬ 
cient transmission between patients in close contact, than is the 
ancestral population from which they emerged. This scenario is 
supported by the observation that hospital-associated MRSA rarely 
establish themselves in the community. Aside from the devel¬ 
opment of antibiotic resistance, little is known about other adap¬ 
tive genetic changes that promote the success of S. aureus (Gomes 
et al. 2005; Chambers and Deleo 2009). A recent microarray-based 
investigation did not reveal any concordant differences between 
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gene contents of highly successful ("epidemic") versus sporadic 
S. aureus strains (Kuhn et al. 2010). Recently, a phage-encoded 
virulence determinant (which is not present in ST22) has been 
shown to be associated with enhanced spreading of MRSA in Asia 
(Li et al. 2012). Phenotypic changes may be caused by genetic 
changes that are subtler than the gain of entire genes through 
horizontal gene transfer, or their loss (Kennedy et al. 2008; Olsen 
et al. 2010). However, our analysis of the genetic changes that ac¬ 
companied the emergence of EMRSA-15 failed to pinpoint a single 
determinant, but, in contrast, identified a series of gene acquisition 
events and mutations that potentially promote the success of this 
clone in the hospital setting (Supplemental Tables S3, S4). Prom¬ 
inent among these are those associated with antibiotic resistance. 

Evidence of the selective pressure exerted by the widespread 
use of broad-spectrum antibiotics in hospitals and beyond is man¬ 
ifest in the ST22 isolates' genomes. In order to unravel the evo¬ 
lutionary events surrounding the development of resistance, we 
examined the genotypic basis of the antibiotic resistance pheno¬ 
types. These analyses permit an evaluation of the utility of genome 
sequencing as a diagnostic tool for the prediction of antibiotic 
resistance. Our in silico analysis of genome sequences identified 
the molecular correlates for >99% of the phenotypic antimicrobial 
resistance traits. This strong correlation documents the impressive 
success of past research into resistance against established antibi¬ 
otics, leaving few undiscovered resistance mechanisms in MRSA. 
Of note, several mutations in the ermC promoter region were 
structurally related, but not identical to previously reported mu¬ 
tations (Levin et al. 2005), and oxazilidinone resistance in our 
isolate collection was due to a mechanism that had been recently 
described from in vitro experiments (Locke et al. 2009), but had 
not previously been reported from clinical isolates (Jones et al. 
2009). These observations highlight the potential of large-scale 
sequence analyses for discovering novel mechanisms. While vari¬ 
ous molecular methods were used in the past to detect resistance 
determinants in clinical isolates (Bergeron and Ouellette 1998), 
whole-genome sequences enabled a screen for all known genes 
and mutations in a single experiment. For the 12 antibiotics, we 
screened for 116 different alleles that were known to confer anti¬ 
biotic resistance by searching antibiotic resistance genes in the de 
novo assemblies or SNPs in the core genome (Supplemental Table 
S5). Accordingly, we were able to demonstrate that it should be 
possible to predict with high accuracy the full resistance profile 
from the sequence data alone. Such predictions may be more 
challenging for other organisms whose mechanisms of resistances 
are poorly studied, or more complex, such as those associated with 
efflux and permeability. 

The reduction in costs and time to generate genomic se¬ 
quences has recently led to genome sequencing being developed as 
a clinical diagnostic tool that can be used to investigate and stop 
outbreaks in hospitals (Eyre et al. 2012; Harris et al. 2012; Koser 
et al. 2012; Snitkin et al. 2012). The observations from this study 
herald an important development in the progress of transitioning 
whole-genome sequencing from a research tool to a diagnostic 
application. Whole-genome sequencing can provide both high 
discriminatory power for molecular epidemiology and a prediction 
of antibiotic resistance, which is valuable for "near real-time" in¬ 
fection control and clinical management, and also for surveillance, 
horizon scanning of the clinical population for the evolution of 
drug resistance, and the identification of emerging threats. 

With the benefit of hindsight it is clear that the emergence of 
the ST22-A2 variant of the EMRSA-15 population in the UK in the 
1980s marked one such threat. Originating in the English Mid¬ 


lands, it spread rapidly during the 1990s, becoming endemic in 
UK hospitals, and was associated with an upsurge in the number of 
MRSA infections (Johnson et al. 2001). 

Of the genetic events that define the emergence of ST22-A2, 
the two nonsynonymous SNPs associated with fluoroquinolone 
resistance would appear to be good candidates to explain the suc¬ 
cess of the clone. Fluoroquinolones are readily excreted in sweat 
and can suppress the growth of the normal microbiota, which 
includes S. aureus that colonizes skin, nose, and throat (Hawkey 
1997). Fluoroquinolone resistant ST22-A2, therefore, may have 
a competitive edge, promoting colonization and survival in the 
hospital environment where this class of antibiotic is frequently 
used. Indeed, widespread fluoroquinolone usage has previously 
been reported to promote the spread of fluoroquinolone-resistant 
MRSA, even though this class of antibiotics is not normally used to 
treat S. aureus infections (Weber et al. 2003; Salangsang et al. 2010). 
A recent study examining MRSA rates in a tertiary hospital over 
a 10-yr period identified a significant increase in MRSA infection 
rates after an increase in fluoroquinolone prescription levels fol¬ 
lowing a period of restricted usage (Parienti et al. 2011). Conversely, 
significant falls in MRSA rates have been linked to a decrease in 
the usage of fluoroquinolones over a 4-yr period (Lafaurie et al. 
2012). In a prevalence study looking at elderly patients admitted 
to hospitals in the UK, a strong association was found between 
MRSA colonization and prior ciprofloxacin exposure (Hori et al. 
2002), with EMRSA-15 being the dominant lineage detected in 
this study (75.9% of MRSA isolates). 

The biological and ecological advantage of fluoroquinolone 
resistance would therefore appear to be clear. However, to under¬ 
stand the remarkable success of ST22-A2 requires insights into the 
historic clinical context in which it emerged. Our analysis suggests 
that development of fluoroquinolone resistance coincided closely 
with the introduction of fluoroquinolone drugs into routine clin¬ 
ical medicine in 1987. Initially, ST22-A2 isolates were restricted to 
the English Midlands until 1990, after which point the clone 
spread rapidly. During this time the use of fluoroquinolones in the 
hospital and community settings increased rapidly, doubling from 
1990 to 1993 (Livermore et al. 2002). Throughout the rest of the 
1990s, usage in the community stabilized, but the usage in hos¬ 
pitals continued to grow steadily both in absolute terms and also as 
a proportion of all antibiotic usage, constituting 31.5% of total 
usage in 1999 compared with 18.9% in 1992 (Livermore et al. 
2002). The emergence and rapid spread of ST22-A2 in the UK in the 
1990s, therefore, perhaps reflects changes in clinical prescription 
regimens. In this regard it is worth noting that we also uncovered 
evidence of the recent influence that country-specific antibiotic 
regimens (clindamycin usage in Germany vs. the UK) have had in 
shaping the MRSA population within Europe, illustrating the link 
between pathogen population structure and national antibiotic 
prescription policymaking. 

The pandemic success of EMRSA-15 cannot, however, be sim¬ 
ply characterized as a case of the right mutation at the right time. 
Fluoroquinolone resistance is not unique to ST22-A2, but is also 
found among non-ST22-A isolates in the collection as the result of 
homoplasic substitutions (Supplemental Table SI), as well as in 
other epidemic MRSA clones (Fluit et al. 2001; Woodford and 
Livermore 2009; Tenover et al. 2011). Hence, other facets of the 
ST22-A2 must have contributed to its outstanding success. Fluo¬ 
roquinolone resistant ST22-A2 emerged from a hospital-adapted 
genetic background (ST22-A1) (Fig. 1) that was already successful, 
and spread through English hospitals by the start of the 1990s 
(Richardson and Reith 1993). 
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An important class of S. aureus proteins involved in host 
colonization are surface-expressed binding proteins (Edwards et al. 
2012). Analysis of the genetic events that define the emergence of 
ST22-A identified mutations in two surface-expressed binding 
proteins, FnBP and Ehb, both of which bind fibronectin, a glyco¬ 
protein component of the host's extracellular matrix. Ebh is one of 
the largest bacterial proteins (—1.1 MDa) (Clarke et al. 2002). The 
pandemic clone of EMRSA-15 has a nonsense mutation that ab¬ 
lates expression of Ehb at the cell envelope. Previous studies dem¬ 
onstrated that mutation of ebh had no effect on the growth rate or 
pathogenicity in a murine skin abscess model of infection (Clarke 
et al. 2002). Notably, this protein is intact in some successful 
MRSA lineages (Gill et al. 2005; Diep et al. 2006; Holden et al. 
2010) and mutated in others (Kuroda et al. 2001; Holden et al. 
2004). Therefore, it is unlikely that mutation of Ehb is the sole 
determinant of increased transmissibility of ST22-A. Its effect may 
be more pronounced due to other mutations, for example, the 
deletion in the fibronectin-binding protein locus. 

Recombination between the C-terminal regions of fhbA and 
fnbB results in an fnbA-fnbB gene fusion and an effective reduction 
in fnb gene copy number at this locus due to the deletion of fnbB. 
Studies with clinical isolates suggested that strains associated with 
invasive disease are significantly more likely to have two fnb genes 
(Peacock et al. 2000). Molecular studies have also shown that both 
FnBPA and FnBPB are required for the establishment of septic in¬ 
fection (Shinji et al. 2011). However, in S. aureus that have both 
homologs, deletion of either gene does not abolish fibronectin- 
binding capacity of the cell (Greene et al. 1995). The functional 
consequences of the deletion and resulting gene fusion on the 
host-cell interaction and transmissibility are unclear, especially 
when considering that there is also a mutation in ebh. 

Interestingly, although most other representatives of the main 
MRSA lineages contain two copies of the fnb genes (Kuroda et al. 
2001; Gill et al. 2005; Diep et al. 2006; Holden et al. 2010), EMRSA- 
16, an unrelated epidemic MRSA lineage that originated in the 
UK around the same time as EMRSA-15, also contains a deletion- 
mediated fibronectin-binding protein gene fusion (Holden et al. 
2004). Like EMRSA-15, EMRSA-16 also has evolved fluoroquinolone 
resistance via gyrA and grlA mutations just prior to its emergence 
(McAdam et al. 2012). Both of these epidemic clones spread 
throughout UK hospitals since the 1990s (Johnson et al. 2001), al¬ 
beit EMRSA-15 has been the more successful of the two clones 
(Ellington et al. 2010). It is therefore intriguing that they share 
the fluoroquinolone resistance and fnbA-fnbB genotypes that 
affect resistance and host cell interaction. However, pinpointing 
the genetic basis of their success is more complex. In the EMRSA- 
16 background, mutations in the accessory gene regulator C 
and a-hemolysin (DeLeo et al. 2011) and squalene desaturase 
(McAdam et al. 2012) have been identified and hypothesized to 
contribute to its epidemic success. However, these are not found 
in EMRSA-15. In order to fully elucidate the genetic basis for the 
relative success of EMRSA-15, extensive experiments investigat¬ 
ing the colonization and transmission properties of the hospital- 
adapted clone are required. 

Methods 

Bacterial isolates 

Sources and properties of 193 S. aureus isolates are listed in Sup¬ 
plemental Table SI. In assembling the collection in this study, we 
included isolates that encompass wide geographical and temporal 


ranges, but also included isolates that reflect the documented 
history and prevalence of EMRSA-15 and ST22. 

Susceptibility to antibiotics was tested by using the stan¬ 
dardized broth microdilution method according to DIN58940 
instructions (Deutsches Institut fuer Normung 2007), and phe¬ 
notypic resistance was defined by applying minimum inhibitory 
concentration breakpoints from the European Committee on An¬ 
timicrobial Susceptibility Testing (EUCAST; http://www.eucast. 
org/). In total, 18 antibiotics were tested (penicillin, PEN; oxacil¬ 
lin, OXA; gentamycin, GEN; linezolid, LNZ; erythromycin, ERY; 
clindamycin, CLI; ciprofloxacin, CIP; fusidic acid, FUS; mupirocin, 
MUP; moxifloxacin, MFL; co-trimoxazole, SXT; tetracyline, TET; 
vancomycin, VAN; teicoplanin,TPL; rifampicin, RAM; phospho- 
mycin, PHO; tigecycline, TGC; daptomycin, DAP), of which 12 had 
resistance in at least one isolate in the collection (Supplemental 
Table SI). Bacterial typing was performed as described previously 
(Strommenger et al. 2008). 

Genomic library creation and multiplex sequencing 

Unique index-tagged libraries for each sample were created, and up 
to 12 separate libraries were sequenced in each of eight channels 
in Illumina Genome Analyzer GAII cells with 54-base paired-end 
reads of Illumina HiSeq with 75-base paired-end reads. The index- 
tag sequence information was used for downstream processing to 
assign reads to the individual samples. 

Detection of SNPs in the core genome 

The paired-end reads were mapped against the chromosome of 
S. aureus HO 5096 0412 (accession no. HE681097) (Koser et al. 
2012), a representative of EMRSA-15 that was isolated from a fatal 
neonatal infection in Suffolk, UK in February 2005. SNPs were 
identified as described in Croucher et al. (2011). Indels were 
identified using Dindel (Albers et al. 2011). Unmapped reads and 
sequences that were not present in all genomes were not con¬ 
sidered as part of the core genome, and therefore SNPs from these 
regions were not included in the analysis. SNPs falling within 
MGEs regions were also excluded from the core genome, as well as 
those falling in high-density SNP regions, which could have arisen 
by recombination. The core genome was curated manually to en¬ 
sure a high-quality data set for subsequent phylogenetic analysis 
and comprised of 2,643,131 bp. 

Phylogenetic analysis 

A phylogeny was drawn for ST22 using RAxML vO.7.4 (Stamatakis 
et al. 2005) to estimate the trees for all SNPs called from the core 
genome. The general time-reversible model with gamma correc¬ 
tion was used for among-site rate variation for 10 initial trees. 

Comparative analysis 

Raw Illumina data were split to generate paired-end reads and 
assembled using a de novo genome assembly program, Velvet 
vO.7.03 (Zerbino and Birney 2008), to generate a multicontig draft 
genome for each of 193 S. aureus isolates. Gene prediction was 
performed using Prodigal (Hyatt et al. 2010). Predicted CDSs were 
compared using BLASTP (Ausubel et al. 1995), and the proteins 
were clustered into homology groups using TribeMCL (Center for 
Mathematics and Computer Science and EMBL-EBI) (Enright et al. 
2002 ) with a cut-off of le -50 . 

Included in the clustering data set were the complete genome 
sequences of S. aureus available in the public sequence databases. 
MGEs in these sequences were identified using comparative 
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genomic analysis and manual curation. Homology groups (HG) 
were categorized according to their origins based on matches to 
CDSs from the MGE reference set that clustered within the group. 
For the purpose of the preliminary accessory genome analysis, HG 
origins were subdivided into Prophage, Conjugative transposon, 
Plasmid, S. aureus Pathogenicity Island (SaPI), SCC mec element, 
Transposon, Hypothetical proteins, and others. Conjugative 
transposon/Plasmid reference sequences included Tn926-like 
(Clewell et al. 1985) and Tn5801 -like elements (Kuroda et al. 2001) 
as well as integrated conjugative elements such as ICE 6013 (Smyth 
and Robinson 2009), whereas Transposon included Tn554 (Phillips 
and Novick 1979), Tn552 (Rowland and Dyke 1989), and IS ele¬ 
ments. Hypothetical protein HGs did not have any matches with 
either core or MGE CDSs from the reference set. HGs defined as 
"Other" indicate that there were no matches with MGE CDSs, but 
there were matches with CDSs from the core regions of the ref¬ 
erence genomes. In some cases, there were examples of clustered 
proteins from more than one origin being present in the single 
HG; for example, Prophage and SaPI and Conjugative trans- 
posons and Plasmid. This reflects the shared evolutionary origins 
and mechanisms of transfer MGEs. In these cases, HGs were as¬ 
signed arbitrarily. 

Detailed comparisons of genome sequences were conducted 
on the de novo assemblies using BLASTN (Ausubel et al. 1995), 
and were facilitated by using the Artemis Comparison Tool (ACT) 
(Carver et al. 2005). 

Determination of antibiotic resistance genotype 

The literature was extensively mined for the genetic basis of anti¬ 
biotic resistance in S. aureus. A database was compiled that con¬ 
tained horizontally acquired genes and core gene mutations that 
conferred antibiotic resistance. Antibiotic resistance conferred by 
SNPs in components of the core chromosome were identified by 
manual curation; substitutions were compared with antibiotic re¬ 
sistance mutations in S. aureus proteins detailed in the literature. 
The presence of antibiotic resistance genes of MGEs was inves¬ 
tigated using BLASTN (Altschul et al. 1990) to compare de novo 
assemblies against a database of S. aureus resistance genes compiled 
from literature and public sequence database searches, and also by 
mapping sequence reads to a pseudo-molecule consisting of con¬ 
catenated antibiotic resistance genes as detailed by Koser et al. 
( 2012 ). 

The genotype of each isolate was searched for the resistance 
genes and/or mutations for the antibiotics tested; for the 12 anti¬ 
biotics for which resistance was detected, 847 incidences in total, 
116 different genes or alleles were searched for. The presence of 
the following known resistance genes and mutations/deletions of 
determinants conferring resistance were identified: mecA (Utsui 
and Yokota 1985), ermC (Catchpole et al. 1988), blaZ (Voladri and 
Kernodle 1998), tetK (Khan and Novick 1983), fusC (O'Neill et al. 
2007), fusA (mutation) (Chen et al. 2010), dfrA (Rouch et al. 1989), 
thy A (deletion) (Leelaporn et al. 1994), ileS-2 (Hodgson et al. 1994), 
aacA-aphD (Rouch et al. 1987), grlA (mutation) (Griggs et al. 2003), 
gyrA (mutation) (Griggs et al. 2003), rplC (mutation) (Locke et al. 
2009), and ermC plasmid (deletion) (Weisblum 1995) (see Supple¬ 
mental Table SI). 

Assessing the association of the number of genetic determinants 
of resistance with time 

The total number of resistance genes and mutations/deletions 
conferring resistance were summed for each strain. Linear re¬ 
gression was used to assess the association of the number of genetic 
determinants of resistance per strain with time. The dependent 


variable was the number of genetic determinants of resistance per 
strain, and the independent variable was time (year of isolation). 
Significance was assessed at P < 0.05. This analysis was applied to 
all 162 ST22-A isolates, from all countries, and separately to the 
77 ST22-A isolates from the UK. This latter analysis was performed 
to evaluate the association of time with the number of genetic 
resistance determinants in a single country. Linear regression with 
the non-ST22-A strains demonstrated that the number of genetic 
determinants of resistance per strain was not significantly associ¬ 
ated with time (P > 0.10). 

Assessing the association of ST22-A strains with hospital origin 

Out of the 162 global ST22-A strains, there were 147 where it was 
known whether or not they were hospital associated: 142 (96%) 
were hospital associated, and five were not. Of the 31 non-ST22-A 
strains, there were 26 where it was known whether or not they 
were hospital associated: 13 (50%) were hospital associated, and 
13 were not. These proportions were significantly different ac¬ 
cording to Fisher's exact test (P < 0.05). 

Bayesian analyses 

We used the Bayesian software package, BEAST (vl.7.1) (Drummond 
et al. 2012) to investigate the temporal, spatial, and demographic 
evolution of ST22. The core genome was retained for the analyses 
after discarding all ambiguous bases and strains that showed ad¬ 
mixture, resulting in a final alignment of 193 taxa and 2,154,246 
bp. The different clades of ST22 (A and non-A) have different de¬ 
mographic histories, as inferred from the shape of their coalescent 
trees. Non-ST22-A individuals coalesce further back in time and at 
a slower rate than ST22-A individuals, and the latter have a more 
star-like genealogy. For this reason, we analyzed these two clades 
separately (25 taxa for non-ST22-A and 162 taxa for ST22-A). We 
used an uncorrelated lognormal distribution (variable clock rate) to 
model the rate of evolution and constrained the tips of the phy- 
logeny to their dates of isolation to calibrate that rate, as well as 
giving it a default uniform prior of between 0 and 1. We used the 
HKY85 model of nucleotide substitution (Hasegawa et al. 1985) 
with four discrete gamma-distributed rate categories, and to model 
the relative node ages, ran the two data sets under both the con¬ 
stant size and exponential size coalescent tree prior. All other de¬ 
fault parameters in the program BEAUti vl.6.2 were used. To model 
the global geographical spread, we calculated the probabilities of 
transfer between each discrete location (i.e., country) via a con¬ 
tinuous-time Markov chain with a nonreversible infinitesimal rate 
matrix, using this to find the location at ancestral nodes (Lemey 
et al. 2009; Edwards et al. 2011). To reduce the number of param¬ 
eters to estimate, we utilized Bayesian stochastic search variable 
selection (BSSVS), where at any point in the Monte Carlo Markov 
chain (MCMC) the transition probabilities of some states can be 
set to zero, with the prior designed to minimize the number of 
location changes as described in Lemey et al. (2009). To infer the 
geographic spread within the UK and to predict the root state of 
the fluoroquinolone-resistant isolates, we modeled diffusion over 
a continuous two-dimensional space assuming a Brownian motion 
(or Weiner process), with the rate of diffusion assumed to be con¬ 
stant over the entire tree and estimated from the data (Lemey et al. 
2010). For this analysis, we used ST22-A strains and discarded 
strains that were fluoroquinolone-sensitive and those that did not 
have exact latitude and longitude information, leaving a total of 
60 taxa. The uncertainty in root location is represented as 80% of 
the highest posterior density of the continuous traits: latitude and 
longitude. We calculated a posterior distribution of all parameters 
using two separate MCMCs and terminated the runs after checking 
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for efficient mixing, convergence, and checking that 90% of the 
distribution was post "burn-in." To estimate the timing of the first 
introduction to each country (Supplemental Fig. SI), we took the 
date of the youngest node for any individual location in the pos¬ 
terior distribution of trees. 

Data access 

The data have been deposited in the European Nucleotide Archive 
(ENA) (http://www.ebi.ac.uk/ena/) under study numbers ERP000150 
and ERP000633. 
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